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The  attenuation  of  a  monochromatic  signal  in  the  presence  of  discrete 
noise  in  one  dimension  is  investigated  numerically.  The  predicted  Gaussian 
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of  a  signal  is  observed.  This  results  shows  that  the  fundamental  property  of  time 
reversibility  can  only  be  established  if  the  overall  system  of  the  waves  and  the 
observer  is  considered.  In  the  second  result,  the  evolution  of  the  amplitude  of  a 
signal  injected  downstream  from  the  noise  is  investigated.  The  Gaussian 
attenuation  is  also  observed  in  this  case.  This  result  explicitly  shows  that  the 
attenuation  length  depends  on  the  distance  the  signal  has  traveled,  thus 
displaying  memory  and  breakdown  of  translational  invariance. 
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I.  INTRODUCTION 


In  the  presence  of  high  intensity  broad  band  noise,  a  small  amplitude 
signal  attenuates  due  to  nonlinear  interactions  with  the  noise.  If  the  noise  is 
isotropic,  the  amplitude  attenuates  as  an  exponential.  This  result  was  first 
predicted  by  Landau  and  Rumer  (1937)  in  connection  with  ultrasonic  attenuation 
due  to  thermal  phonons  in  a  dielectric  at  low  temperature.  Westervelt  (1976) 
extended  these  ideas  to  the  case  where  the  broad  band  noise  was  due  to 
classical  sound  waves.  Experimental  observations  of  the  attenuation  of  sound  in 
superfluid  helium  (Maris,  1973)  and  of  sound  in  water  (Stanton  and  Beyer,  1978 
and  1981)  support  these  results. 

The  limiting  case  of  interaction  with  noise  in  one  dimension  is  derived 
from  Riemann’s  exact  equation  for  lossless  propagation  of  unidirectional  sound 

Ut+(c+pu)u^  =  0  ,  (1.1) 

where  u  is  the  particle  velocity,  p  =  1  +  (p/c)(dc/dp)  is  the  nonlinear  coefficient,  c 
is  the  equilibrium  value  of  the  speed  of  sound,  and  p  is  the  ambient  density.  The 
amplitude  of  the  signal  is  calculated  to  attenuate  as  a  Gaussian  exp(-rx^),  where 
X  is  the  distance  from  the  source,  and  where  the  nonlinear  attenuation  coefficient 
is 


t 


(1.2) 


where  fo  is  the  frequency  of  the  signal  and  If  is  the  spectral  intensity  of  the  noise 
at  frequency  f  (Rudenko  and  Chirkin,  1975;  Rudenko  and  Soluyan,  1977).  The 


1 


integral  is  the  intensity  of  the  noise.  Hence,  in  the  shockless  regime  3  is 
independent  of  the  details  of  noise  spectral  distribution.  It  should  be  noted  that 
this  is  consistent  with  the  exp(-rx^)  solution  to  the  evolution  equation  dA/dx  =  - 
2rxA  for  the  amplitude  A.  If  r  were  dependent  upon  the  distribution,  then  the 
spectral  evolution  of  the  noise  due  to  nonlinear  interactions  would  necessarily 
cause  r  to  be  a  function  of  position,  which  would  invalidate  the  Gaussian 
solution. 

The  basic  nonlinear  interaction  of  acoustic  waves  is  a  three-wave 
resonance,  so  that  waves  with  frequency  fi  and  fa  scatter  to  produce  waves  with 
frequencies  fa  =  |  h  ±  fa  I .  Furthermore,  because  these  waves  are  nondispersive 
the  interactions  can  only  be  collinear  (Rudenko  and  Soluyan,  1977).  For 
random  waves  in  two  and  three  dimensions,  energy  in  the  noise  components 
redistributes  irreversibly  by  collinear  resonant  interactions.  Randomization  in 
angle  is  brought  about  by  slower  collision  processes  with  a  local  transfer 
between  adjacent  rays  (Newell  and  Aucoin,  1971).  After  a  time  determined  by 
the  nonlinearities,  a  closed  system  of  interacting  random  waves  can  thus  reach 
thermodynamic  equilibrium,  and  small  fluctuations  about  this  state  relax 
exponentially  with  an  attenuation  coefficient  proportional  to  the  noise  intensity. 
Because  the  noise  effectively  acts  as  a  heat  bath,  it  is  appropriate  to  refer  to  the 
absorption  of  sound  by  noise  in  this  case. 

The  predicted  Gaussian  attenuation  in  one-dimensional  nonlinear 
acoustics  brings  together  two  physically  significant  results.  First,  because  a 
small  amplitude  signal  can  be  regarded  as  a  fluctuation  of  the  noise,  and 
because  this  fluctuation  does  not  relax  exponentially,  we  infer  that  a  system  of 
interacting,  collinear,  nondispersive  waves  is  inherently  far  off  equilibrium.  This 
results  from  both  a  lack  of  angle  randomization  and  the  fact  that  all  spectral 
components  are  in  resonance  with  each  other  in  this  case. 
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A  second  aspect  of  a  Gaussian  attenuation  is  that  it  breaks  translational 
invariance  because  the  attenuation  coefficient  is  now  a  function  of  position. 
Specifically,  by  taking  measurements  at  several  distances  from  the  source  of  a 
signal,  one  can  calculate  both  the  location  and  strength  of  the  source.  This 
cannot  be  done  if  the  attenuation  is  exponential.  This  is  physically  significant 
because  the  equations  that  describe  the  evolution  of  nonlinear  noise  and  its 
interaction  with  a  monofrequency  signal  are  translationally  invariant. 

An  experiment  to  measure  the  predicted  Gaussian  attenuation  has  been 
carefully  investigated  by  Larraza  et  al.  (1996).  In  this  experiment,  the 
dimensionality  is  controlled  by  performing  measurements  in  a  long  traveling 
wave  tube  driven  at  frequencies  below  the  first  cutoff.  The  theory  is  in  excellent 
agreement  with  the  experimental  results  only  after  it  is  modified  to  incorporate 
wall  losses.  The  agreement  is  shown  as  the  frequency,  noise  level,  and 
distance  from  the  source  are  varied.  In  addition,  Larraza  et  al.  (1996)  observed 
the  spectral  intensity  of  the  high-frequency  tail  of  fully  developed  shockless 
noise  to  be  an  f"^  power  law  in  the  frequency  f,  in  accord  with  Kutznetsov’s 
theory  (1970).  This  power  law  is  a  consequence  of  the  far  off  equilibrium  nature 
of  the  system. 

The  purpose  of  this  thesis  is  to  numerically  investigate  the  predicted 
Gaussian  attenuation.  In  Ch.  II  we  give  a  general  theoretical  background  for 
nonlinear  acoustics  in  one  dimension.  In  Ch.  Ill  we  describe  the  details  of  the 
numerical  code.  The  program  is  tested  by  comparing  the  numerical  results  with 
the  predictions  of  theory  for  the  cases  of  pure  tone  propagation  and  suppression 
of  a  signal  by  a  pump.  In  Ch.  IV  we  present  numerical  simulations  of  the 
absorption  of  sound  by  noise.  We  verify  the  absorption  of  sound  by  noise  and 
also  present  two  new  results.  In  the  first,  we  investigate  the  transition  from 
Gaussian  to  Bessel  as  a  function  of  resolution  in  the  detection  of  a  signal  in  the 
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presence  of  noise.  In  the  second  result,  as  a  test  of  the  breakdown  of 
translational  invariance,  we  investigate  the  evolution  of  the  amplitude  of  a  signal 
injected  downstream  in  the  presence  of  noise.  Conclusions  and  possible  future 
work  are  presented  in  Ch.  V. 
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il  THEORY 


In  this  chapter,  we  describe  Riemann’s  equation  for  sound  propagation  in 
one  dimension  in  a  barotropic  fluid,  and  then  derive  the  general  shock  inception 
distance.  We  also  present  the  analytic  solutions  to  some  special  cases.  The 
simplicity  of  Riemann’s  solution  allows  the  determination  of  exact  and  asymptotic 
explicit  solutions  in  the  preshock  region.  Examples  include  Fubini’s  (1935)  pure 
tone  radiation  problem,  Fenlon’s  (1970)  generalization  to  multiple  frequency 
sinusoidal  source,  Kutznetsov’s  (1970)  asymptotic  f“®  spectrum  for  broadband 
noise,  and  Rudenko  and  Chirkin’s  (1975)  Gaussian  attenuation  of  a  signal  due  to 
interactions  with  normally  distributed  random  noise. 


A.  RIEMANN’S  EQUATION 


The  basic  equations  of  motion  for  an  isentropic  perfect  gas  are  governed 
by  the  continuity  equation  and  Newton’s  second  law 

Pt  +  ^Px  +  P^x  =  0  >  (2-1a) 

Ut  +  uux  +  p'Vx  =  0  .  (2.1b) 

respectively,  where  p  is  density,  p  is  total  pressure,  and  u  is  the  particle  velocity. 
The  subscripts  x  and  t  denote  partial  differentiation  with  respect  to  the  space  and 
time  coordinates,  respectively.  Eqs.  (2.1 )  must  be  supplemented  with  the 
equation  of  state  of  an  isentropic  ideal  gas 


Po 


(2.2) 
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where  po  and  Po  are  the  ambient  density  and  pressure,  respectively,  and  y  is  the 
ratio  of  specific  heats  of  the  gas. 

For  the  situation  where  the  waves  are  propagating  in  one  direction,  Eqs. 
(2.1)  and  (2.2)  imply  Riemann’s  equation  (Blackstock,  1972) 

Ut +(c  +  pu)u^  =0  ,  (2.3) 

where  c  =  {dpldp\  is  the  equilibrium  value  of  the  speed  of  sound,  and  where  for 
an  ideal  gas,  the  nonlinear  coefficient  is  given  by  p  =  (1+y)/2.  From  Eq.  (2.3),  the 
instantaneous  speed  of  propagation  of  a  point  in  the  wave  with  particle  velocity  u 
is  c+pu.  Thus,  points  of  the  wave  having  different  particle  velocity  propagate  at 
different  speeds. 

A  general  implicit  solution  of  Eq.  (2.4)  with  the  x  =  0  boundary  condition 
u(0,t)  =  f(t)  is  given  by 


u(x,t)  =  f 


^  X  ^ 

t - 

V  c  +  pu(x,t)> 


(2.4) 


as  can  be  verified  by  direct  substitution  into  Eq  (2.3).  Thus  the  velocity  at  the 
boundary  x  =  0  determines  the  velocity  at  x  >  0  through  the  retarded  time 


T  =  t  - 


X 

C  +  pu  ■ 


(2.5) 


B.  SHOCK  INCEPTION  DISTANCE 

Suppose  Ui  is  produced  at  x  =  0  at  time  t^,  and  a  greater  velocity  U2  is 
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produced  at  a  later  time  t2.  Riemann’s  equation  breaks  down  when  the  second 
disturbance  overtakes  the  first,  which  will  eventually  occur  at  some  x  >  0.  This 
overtaking  yields  a  discontinuity  (infinite  slope  or  shock)  in  u.  If  t^  and  ts  differ  by 
a  finite  amount,  a  discontinuity  will  occur  at  a  distance  less  than  x,  because  the 
waveform  u(x,t)  continuously  distorts  as  time  progresses.  We  are  thus  led  to 
consider  velocity  disturbances  u  and  u+du  and  at  x  =  0  that  are  an  infinitesimal 
time  dt  apart.  If  du  >  0,  the  disturbances  will  coincide  at  some  distance  L.  If  T  is 
the  elapsed  time. 


[c  +  pu]T  =  L  , 

(2.6a) 

[c  +  p(u  +  du)](T  -  dt)  =  L  , 

(2.6b) 

corresponding  to  the  first  and  second  disturbances,  respectively.  Multiplying  out 
Eq.  (2.6b),  substituting  Eq.  (2.6a),  and  solving  for  L,  yields 


(c  +  pf) 

pdf/dt  ■ 


The  quantities  f  and  df/dt  are  the  values  of  the  velocity  and  its  slope  at  x  =  0, 
respectively.  The  expression  (2.7)  gives  the  distance  at  which  an  initial  positive 
slope  becomes  infinite.  (A  negative  slope  never  becomes  infinite.)  For  any 
solution  to  Riemann’s  equation  to  be  valid  at  some  distance  x,  it  must  be  verified 
that  X  <  L  for  all  values  of  f(t). 


C.  PURE  TONE  PROPAGATION 

For  the  propagation  of  a  pure  tone,  the  solution  (2.4)  subject  to  the  x  =  0 
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boundary  condition 


u(0,t)  =  Uq  sin(a)t)  , 

may  be  written  in  the  form  of  Fourier  series 


(2.8) 


u  =  Uo^BnSin(n(Bt)  , 
n=1 


(2.9) 


where  in  the  preshock  region  the  Fourier  coefficients  are  given  by  (Fubini,  1935) 


2Jr,(pUocox/c^) 

npuocox/c^ 


(2.10) 


where  is  the  Bessel  function  of  order  n. 

For  the  x  =  0  velocity  (2.8),  the  minimum  value  of  the  shock  length  (2.7) 
corresponds  to  cot  being  an  integral  multiple  of  27t.  The  minimum  distance  at 
which  shocking  occurs  for  a  pure  tone  is  thus 


L  = 


PUqCO 


(2.11) 


D.  SUPPRESSION  OF  A  SIGNAL  BY  A  PUMP 

The  simplest  problem  of  interaction  of  sound  with  sound  is  given  by  the  x 
=  0  boundary  condition 
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u(0,t)  =  Up  sin^copt)  +  Ug  sin(o)st)  , 


(2.12) 


where  we  focus  here  on  the  interaction  of  a  strong  pump  wave  of  low  frequency 
oip  with  a  weak  signal  of  high  frequency  cOg.  The  signal  has  little  effect  on  the 
propagation  of  the  pump.  However,  the  pump  will  modulate  the  signal  and 
generate  sidebands  at  frequencies  cogtcop,  thereby  removing  energy  from  the 
signal.  An  approximation  to  Fenlon’s  (1970)  solution  for  the  evolution  of  the 
signal  is 


u(x,t)  =  UsJo(pUpC0sx/c2)sin[®s(t-x/c)]  .  (2.13) 

The  amplitude  of  the  signal  depends  only  on  the  particle  velocity  of  the  pump, 
the  frequency  of  the  signal,  and  the  distance.  The  signal  vanishes  at  the  zeros 
of  the  Bessel  function,  at  which  the  energy  has  been  temporarily  completely 
pumped  into  adjacent  sidebands. 


E.  ABSORPTION  BY  NOISE 


For  the  case  of  a  weak  signal  of  frequency  f  in  the  presence  of  finite- 
amplitude  broadband  noise,  the  amplitude  of  the  signal  attenuates  as  a 
Gaussian  exp(-rx^),  where 


r  = 


27t2p2f2 

- 

4  '^rms  > 


(2.14) 


where  u^^g  is  the  rms  velocity  of  the  noise  (Rudenko  and  Chirkin,  1975;  Rudenko 
and  Soluyan,  1977).  Broadband  noise  may  be  considered  by  the  x  =  0  boundary 
condition 


9 


N 

^noise(®'^)  ~  9n)  >  (2-15) 

n=1 

where  A  is  the  peak  velocity  of  each  component  of  the  noise,  (o^  are  densely 
distributed  frequencies,  and  cpn  are  randomly  distributed  phases.  The  Gaussian 
attenuation  thus  results  as  a  limiting  case  of  multiple  pump  waves  with  random 
phases.  When  a  signal  is  interacting  with  one  pump,  the  restitution  effect 
(amplitude  increasing  with  distance)  is  lost  when  the  interaction  occurs  with  a 
large  number  of  pumps  with  random  phases.  Because  of  the  central  limit 
theorem,  the  noise  described  by  (2.15)  should  approach  a  normal  distributed 
noise  for  a  sufficiently  large  number  N.  In  this  limit,  a  theory  formulated  in  terms 
of  the  noise  (2.15)  should  accurately  yield  the  Gaussian  attenuation. 

By  squaring  Eq.  (2.15)  and  averaging  overtime,  we  determine  the  square 
rms  particle  velocity  of  the  noise.  Solving  for  A  yields 

(2.16) 


The  approach  to  the  Gaussian  attenuation  can  be  probed  by  increasing  N  while 
keeping  Ur^s  constant. 
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III.  NUMERICAL  METHOD 


In  this  chapter,  we  describe  the  numerical  method  used  to  integrate 
Riemann’s  equation  for  shockless  sound  propagation  in  one  dimension.  We  test 
the  numerical  program  by  comparing  its  results  to  theory  for  the  cases  of  pure 
tone  radiation  and  suppression  of  a  signal  by  a  pump. 

A.  FORWARD  PROPAGATION  METHOD 

As  a  numerical  convenience,  we  choose  p  =  1  and  c  =  1  so  that  the 
dimensionless  Riemann’s  equation  (2.3)  and  the  expression  (2.5)  for  the 
retarded  time  become 

du  \  5u  „ 

■^+(l+u)-  =  0,  (3.1) 


and 


(3.2) 


respectively.  As  stated  in  Sec.  2.1,  a  general  implicit  solution  of  Eq.  (3.1)  with 
the  X  =  0  boundary  condition  u(0,t)  =  f(t)  is  given  by  u(x,t)  =  f(T). 

Eq.  (3.2)  determines  a  set  of  straight  lines  in  the  x-x  plane.  The  slope  of 
each  line  is  determined  by  the  value  of  the  velocity  at  the  boundary  at  a  given 
time  (wavelet).  These  lines  are  known  as  characteristics,  where  the  parameter  t 
is  the  time  base  for  the  source  signal  with  t  =  x  at  x  =  0. 

The  principle  of  the  forward  propagation  method  stems  from  the  fact  that 
the  characteristics  represent  the  trajectories  of  motion  of  the  different  wavelets  in 
the  x-x  plane.  The  velocity  of  propagation  of  a  wavelet  is  equal  to  1+u  where  u 
is  the  assigned  value  of  the  velocity  of  the  wavelet  at  the  boundary.  The  method 
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is  implemented  by  storing  at  the  boundary  the  wavelets  and  their  corresponding 
times.  In  the  propagation  of  the  initial  waveform  the  first  incremental  distance 
produces  a  new  time  array  associated  with  the  velocity  array.  While  at  the  origin 
the  spacing  between  time  arrays  is  commensurate,  once  the  waves  propagate  a 
small  distance  the  time  between  adjacent  elements  is  stretched  or  compressed 
according  to  Eq.  (3.2),  which  physically  corresponds  to  distortions  of  the  initial 
waveform.  Appendix  A  shows  a  C-language  numerical  program  that  implements 
the  forward  propagation  method. 

B.  DISCRETE  FOURIER  TRANSFORM 

Our  main  emphasis  is  in  determining  the  evolution  of  the  amplitude  of  a 
monochromatic  signal  in  the  presence  of  noise.  As  we  forward  propagate  the 
boundary  condition  of  a  time  series  with  the  signal  plus  noise,  sidebands  of  the 
noise  appear  about  the  signal’s  frequency.  Also,  the  amplitude  of  the  signal  is 
expected  to  attenuate  with  distance.  It  is  thus  necessary  to  extract  the  signal 
from  the  noise  at  a  distance  x  from  the  source.  We  achieve  this  by  Fourier 
analyzing  the  waveform  at  the  frequency  of  the  signal,  which  amounts  to  time 
averaging  the  in-phase  and  quadrature  components  of  the  waveform.  This 
process  is  referred  to  as  the  “discrete  Fourier  transform”  (DFT)  because  the 
integration  is  performed  with  a  finite  (rather  than  continuous)  time  increment  and 
is  performed  over  a  finite  (rather  than  infinite)  time  interval. 

The  integration  in  our  DFT  is  accomplished  with  trapezoidal  summation. 
(Simpson’s  rule  cannot  be  employed  because  the  time  increment  is  not 
constant.)  Appendix  A  shows  the  numerical  program  which  includes  the 
trapezoidal  summation. 

The  purpose  of  this  section  is  to  examine  the  effect  of  a  finite  rather  than 
infinite  time  interval.  In  the  first  part,  we  examine  the  error  that  arises  when  a 
time  series  is  averaged  over  a  finite  interval.  We  will  then  apply  the  results  to 
the  case  of  determining  the  amplitude  of  a  harmonic  signal  in  the  presence  of 
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equally-spaced  sidebands.  Proper  choice  of  the  averaging  time  is  shown  to  lead 
to  the  complete  elimination  of  error. 

1.  Averaging  Over  Finite  Time  intervals 

Consider  any  time  series  F(t).  This  can  be  represented  by  the  Fourier 
expansion 


F(t)  =  ^[G(f)cos(27rft)  +  H(f)sin(27tft)]  ,  (3.3) 

f 

where  f  is  the  frequency.  The  dc  component  of  F(t)  can  be  analytically  obtained 
by  averaging  over  an  infinite  time: 

T 

G(0)  =  ^  jF(t)dt.  (3.4) 

-T 

For  finite  values  of  the  time  T,  the  average  only  yields  an  approximate 
value  of  the  dc  component.  The  general  finite-time  average  is 


F  = 


(3.5) 


where  To  is  arbitrary.  Substituting  the  Fourier  expansion  (3.3)  into  the  average 

(3.5),  and  performing  the  integration,  yields 


F 


sin[27tf(T  +  To )]  -  sin(27tfTo ) 
2^cff 


G(f) 


(3.6) 


cos[27if(T  +  To )]  -  cos(27ifTo ) 
27rfT 


In  the  limit  T  ->  oo,  note  that  Eq.  (3.6)  correctly  reduces  to  G(0).  Eq.  (3.6)  shows 
that  the  effect  of  averaging  over  a  finite  time  interval  is  to  low-pass  filter  the  time 
series.  The  filtering  depends  upon  the  phase  of  the  time  window  compared  to 
the  time  series.  For  example,  for  an  average  over  -T/2  to  112  (i.e.,  To  =  -T/2), 
the  filter  applied  to  G(f)  is  sin(7ifT)/(7ifT),  and  the  filter  applied  to  H(f)  is  zero. 

This  is  in  fact  the  simplest  case  with  which  to  deal,  but  will  continue  with  our 
general  choice  (3.5)  in  order  to  ensure  the  generality  of  the  results. 

In  Eq.  (3.6),  note  that  the  filters  simultaneously  vanish  if  27rfT  =  2%,  4%,  6k, 
... ,  except  for  the  G(f)  filter  at  f  =  0.  Hence,  if  the  Fourier  components  of  the  time 
series  are  restricted  to  the  equally-spaced  values  f  =  1/T,  2/T,  3/T,  ... ,  no  error 
results  in  the  finite-time  average]  i.e.,  the  value  of  Eq.  (3.5)  is  G(0).  The  reason 
for  this  is  simply  that  the  integration  interval  T  in  Eq.  (3.6)  corresponds  to 
integral  numbers  of  cycles  of  these  frequencies,  so  they  will  all  average  to  zero. 

2.  Determination  of  the  Amplitude  of  a  Signal  in  the  Presence  of 
Discrete  Noise 

We  now  consider  the  time  series 

F(t)  =  Acos(27tft  + y)  +  Bcos[27t(f  ±  Af)t  +  (p]  .  (3.7) 

The  first  term  is  considered  to  be  the  “signal,”  and  the  second  term  the  “noise,” 
where  Af  >  0.  The  problem  is  to  determine  the  signal  amplitude  A  from  F(t) 
alone.  Analytically,  this  can  be  accomplished  by  down-converting  the  signal 
frequency  f  to  dc,  and  then  averaging; 
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(3.8a) 


T 

Ac  =  lim  ^  fF(t)cos(27tft)dt  , 

T->00  T  J 
-T 

T 

As  =  lim  ^  fF(t)sin(27tft)dt  ,  (3.8b) 

T^oo  I  J 
-T 

A  =  (a,  +  A,y“  .  (3.8c) 

Compared  to  the  dc  case  (3.4),  a  factor  of  2  has  been  included  in  Eqs.  (3.8a, b) 
because  half  of  the  amplitude  at  frequency  f  is  converted  to  zero  frequency. 

(The  other  half  is  converted  to  frequency  2f.) 

As  examined  in  Sec.  3.2A,  the  integrations  (3.8a,b)  will  in  general  only  be 
approximate  for  a  finite  averaging  time  T.  We  thus  consider  the  now- 
approximate  quantities 


Ac  =  Y  |F(t)cos(27tft)dt  ,  (3.9a) 

T„ 


As  =  Y  J'F(t)sin(27ift)dt  ,  (3.9b) 

To 

where  To  is  arbitrary. 

The  first  problem  is  to  ensure  that  the  results  are  exact  when  B  =  0  in  Eq. 
(3.7).  This  is  readily  seen  to  occur  if  the  averaging  time  is  an  integral  number  of 
cycles  of  the  signal; 
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T 


n 

T 


n  =  1,2,  3,  ...  . 


(3.10) 


As  shown  in  Sec.  3.2A,  averaging  over  a  finite  time  interval  has  the  effect 
of  introducing  a  low-pass  filter  that  yields  no  error  if  the  frequency  spacing  of  the 
unwanted  Fourier  components  is  1/T  from  zero  frequency.  Due  to  the  down¬ 
converting  in  Eqs.  (3.9a, b),  this  corresponds  here  to  unwanted  Fourier 
components  with  spacing  1/T  about  the  frequency  f  of  the  signal.  That  is,  there 
is  no  error  due  to  a  component  of  frequency  f  ±  mAf,  where  m  =  1 ,  2,  3,  . . . ,  if 


Combining  Eqs.  (3.10)  and  (3.1 1 ),  we  find  that  all  error  will  be  eliminated  if 
f 

n  =  ^.  n  =  1,2, 3,  ...  (3.12) 

In  the  implementation  of  a  DFT,  the  integration  time  n/f  and  the  frequencies  f  and 
Af  must  be  chosen  such  that  the  condition  (3.12)  holds.  In  practice,  the 
frequencies  are  usually  chosen  to  be  convenient  numbers,  so  that  their  ratio  is 
automatically  an  integer.  The  result  (3.12)  is  explicitly  verified  in  the  next 
section. 

3.  Explicit  Calculation  of  the  Error  Due  to  Finite-Time  Averaging  of 
the  Noise 

From  Eqs.  (3.9a,b),  the  errors  that  arise  due  to  the  presence  of  the  noise 
in  Eq.  (3.7)  are 
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AA,  =  Y 


I  T  Iq 

j*  cos[27i(f  ±  Af)t  +  (p]cos(27tft)dt 


(3.13a) 


2B  °  ° 

—<  c»s((p)  I  cos[2ii(f  ±  Af)t]cos(27ift)dt  -  sin((p)  J  sin[27c(f  ±  Af)t]cos(27tft)dtl  . 


AAs  =  Y 


I  -r  Iq 

J  cos[27c(f  ±  Af)t  +  (p]sin(27ift)dt  , 


(3.13b) 


I  +  lo  '  +  'o 

Y'  cos((p)  I cos[27i(f  ±  Af)t] sin(27tft) dt  -  sin((p)  | sin[27t(f  ±  Af)t]sin(2i:ft)dt 


Using  the  identities  cos(x)cos(y)  =  [cos(x+y)+cos(x-y)]/2,  sin(x)cos(y)  = 
[sin(x+y)+sin(x-y)]/2,  and  sin(x)sin(y)  =  [-cos(x+y)+cos(x-y)]/2,  we  express  Eqs. 
(3.13a,b)  as 

t+t, 

AAj.  =  —cos{<p)  J  {cos[27r(2f  ±  Af)t]  +  cx)s(27cAft)}  dt 


(3.14a) 


B  ° 

-  YSin((p)  J {sin[27i(2f  ±  Af)t]  +  sin(27cAft)}  dt 


AA,  =  - 


'  *0 

cos((p)  J  {sin[27i(2f  ±  Af)t]  -  sin(27tAft)}  dt 


(3.14b) 
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+ 


B 

T 


T+T„ 

sin((p)  J  {cos[27i(2f  ±  Af)t]  -  cos(27cAft))  dt  . 

To 


Consider  the  second  terms  in  the  integrands  in  Eqs.  (3.14a,b).  Due  to  Eq. 

(3.1 1 ),  the  integration  is  over  one  cycle  of  these  terms,  so  the  integration  yields 
zero.  Using  Eqs.  (3.10)  and  (3.1 1 ),  we  can  express  the  period  of  the  first  terms 
as  1/(2f  ±  Af)  =  1/(2n/T  ±  1/T)  =  T/(2n  ±  1 ).  Because  this  integration  time  T  is  an 
integral  number  of  this  period,  the  integration  yields  zero.  Hence,  the  errors 
vanish  (AAc  =  AAs  =  0),  and  so  sidebands  at  f  ±  Af  do  not  contribute  to  Eqs. 
(3.9a,b). 

Now  suppose  that  there  is  another  sideband  at  f  ±  2Af,  and  that  Eqs. 
(3.10)  and  (3.1 1)  still  hold.  The  second  terms  in  the  integrands  of  Eqs.  (3.14a,b) 
still  vanish,  because  the  integrations  are  over  two  cycles.  The  period  of  the  first 
terms  is  now  T/[2(n  ±  1 )],  which  still  leads  to  vanishing  integrals  because  the 
integration  time  T  is  an  integral  number  of  this  period.  The  argument 
generalizes  to  sidebands  at  f  ±  mAf,  where  m  =  1 ,  2,  3 . 

C.  TESTS  OF  THE  NUMERICAL  PROGRAM 


To  test  the  numerical  program,  we  consider  two  simple  cases.  One  is  the 
propagation  of  a  pure  tone  and  the  other  is  the  suppression  of  a  weak  signal  of 
high  frequency  by  a  high  amplitude  pump  of  low  frequency  in  the  preshock 
region.  In  both  cases,  the  numerical  results  are  in  excellent  agreement  with 
analytical  theories. 

1.  Pure  Tone  Propagation 

Consider  the  pure  tone 
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u(0,t)  =  Uq  sin(27cft) 


(3.15) 


at  X  =  0.  Due  to  nonlinear  effects,  energy  flows  to  harmonics  of  the  fundamental 
frequency  which  therefore  grow  as  a  function  of  distance.  In  the  preshock 
region,  the  energy  in  the  higher  harmonics  equals  the  energy  loss  at  the  main 
frequency.  Consequently,  the  amplitude  at  frequency  f  decreases  with  distance. 
From  Eq.  (2.10),  in  dimensionless  variables  (c  =  p  =  1)  the  amplitude  as  a 
function  of  distance  is 


2Ji(27CfUoX) 

27rfUoX 


(3.16) 


According  to  the  shock  distance  (2.1 1 ),  Eq.  (3.16)  is  valid  only  for  distances  x  < 
1/27rfuo,  and  the  amplitude  thus  decreases  monotonically. 

Figure  3.1  shows  a  comparison  between  the  numerical  integration  and 
theory,  for  a  signal  with  amplitude  Uo  =  10~'^  and  frequency  f  =  50.  The  shock 
inception  distance  equals  31 .83  units  in  this  case.  For  20  time  steps  per  cycle 
there  is  a  small  but  significant  error.  Increasing  the  number  of  time  steps  to  200 
gives  excellent  agreement  with  theory. 
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Fig  3.1  Normalized  amplitude  as  a  function  of  distance  for  a  pure  tone  injected  at  the  origin. 

The  amplitude  of  the  pure  tone  is  1 0  and  the  frequency  is  50.  The  points  are  the  result 
of  numerical  integration,  and  the  solid  line  corresponds  to  the  theoretical  expression 
(3.16).  The  squares  correspond  to  20  time  steps  per  cycle,  and  the  circles  correspond  to 
200  time  steps  per  cycle. 


2.  Suppression  of  a  Signal  by  a  Pump 

As  discussed  in  Ch.  II.  and  Sec.  D.,  the  suppression  of  a  signal  by  a  pump 
is  formulated  in  terms  of  the  boundary  condition  (2.12)  at  x  =  0 

u(0.  t)  =  Up  sin(27rfpt)  +  Ug  sin(27tf3t)  ,  (3.17) 

where  the  pump  frequency  is  fp  =  1  and  the  signal  frequency  is  fs  =  50.  The 
weak  signal  has  amplitude  10"'*  and  has  little  effect  on  the  propagation  of  the 
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pump  with  amplitude  10"^.  However,  the  pump  will  modulate  the  signal  and 
generate  sidebands  at  frequencies  +  fp,  thereby  removing  energy  from  the 
signal.  In  dimensionless  variables,  Fenlon’s  approximate  solution  (2.13)  for  the 
evolution  of  the  signal  is 


UsJo 


(27tfsUpx) sin[27tfs 0  -  x  /  c)]  . 


(3.18) 


Figure  3.2  shows  a  comparison  between  the  numerical  integration  and  the 
theory  (3. 1 8).  The  shock  inception  distance  (2.1 1 )  for  the  pump  equals  1 5.92 
units.  The  number  of  time  steps  per  cycle  of  the  signal  is  20.  There  is  excellent 
agreement  between  theory  and  numerical  integration. 


Fig  3.2  Suppression  of  a  signal  with  frequency  50  by  a  pump  of  amplitude  10“^  and  frequency  1. 
The  signal’s  amplitude  is  normalized  to  its  value  of  10"'*  at  the  origin.  The  points  are  the 
result  of  numerical  integration,  and  the  solid  line  corresponds  to  the  theoretical 
expression  (3.18).  The  number  of  time  steps  per  cycle  of  the  signal  is  20. 
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IV.  NUMERICAL  SIMULATIONS  OF  ABSORPTION  OF  SOUND  BY 

NOISE 


In  this  chapter  we  present  numerical  simulations  of  the  absorption  of 
sound  by  noise  in  one  dimension.  We  verify  the  Gaussian  attenuation  of  a 
signal  in  the  presence  of  noise  and  present  two  new  results.  In  the  first,  we 
investigate  the  transition  from  Gaussian  to  Bessel  dependence  as  a  function  of 
resolution  in  the  detection  of  a  signal.  In  the  second  result,  as  a  test  of  the 
breakdown  of  translational  invariance,  we  investigate  the  evolution  of  the 
amplitude  of  a  signal  injected  downstream  from  the  noise. 

A.  ABSORPTION  BY  NOISE 

As  stated  in  Ch.  I.  and  Sec.  E.,  a  small-amplitude  signal  of  frequency  f 
interacting  with  finite-amplitude  broadband  noise  is  predicted  to  attenuate  in 
distance  as  a  Gaussian.  We  choose  the  signal  to  be  injected  with  peak  particle 
velocity  Us  at  x  =  0,  so  the  particle  velocity  of  the  signal  at  x  =  0  is 

Usignal(O.t)  =  U3COS(27tft)  .  (4.1) 

The  rms  (root-mean-square)  amplitude  of  the  particle  velocity  of  the  signal  as  a 
function  of  position  is  predicted  to  attenuate  as 

Vsignai  =  ;^exp(-rx2)  ,  (4.2) 

where  the  Gaussian  coefficient  (2.14)  is,  in  dimensionless  units, 

r  =  2n^f^L,  (4.3) 
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where  Urms  is  the  dimensionless  rms  particle  velocity  of  the  noise,  which  is 
constant  in  space  (although  the  spectral  components  of  noise  change). 
The  noise  is  described  at  x  =  0  by  the  particle  velocity 


*^noise(®i 


Ur 


Z  cos 


2% 


n-1 


1+  M 

N 


(4.4) 


where  cpn  are  random  phases.  Note  that  the  total  average  acoustic  intensity  of 
the  noise  (4.4)  is 


u 


2  _ 
rms  ” 


(4.5) 


The  frequency  components  of  the  noise  (4.4)  are  equally  spaced  in  a 
band  of  frequencies  between  1  and  2.  In  the  preshock  region,  the  theory  should 
be  accurate  (i.e.,  the  signal  should  attenuate  as  the  predicted  Gaussian)  for  a 
sufficiently  large  number  N  of  components  of  the  noise.  If  all  the  phases  in  Eq. 
(4.4)  are  set  to  kI2,  the  noise  shock  inception  distance  is  (Fenlon,  1970) 


L  = 


Vn 

7cUo(3N-1)  ■ 


(4.6) 


For  amplitude  uo  =  0.01  and  N  =  50  components,  the  shock  inception  distance 
(4.6)  is  approximately  L  =  1.51  units.  According  to  Eq.  (2.7)  random  phases 
would  yield  a  larger  shock  inception  distance  because  Eq.  (4.4)  with  all  the 
phases  set  to  n/2  has  the  maximum  slope  at  the  origin. 

Figure  4.1  shows  a  comparison  between  the  numerical  integration  and 
theory,  for  a  signal  with  amplitude  Us  =  1 0  ''  and  frequency  f  =  50.  In  the 
numerical  integration,  there  are  20  time  steps  per  cycle  of  the  signal.  Note  that 
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the  amplitude  of  the  signal  appears  to  asymptote  to  a  nonzero  constant  value. 
Increasing  the  averaging  time  does  not  improve  this  condition.  Also,  for  the 
parameters  of  the  problem,  we  numerically  determined  the  shock  inception 
distance  to  be  greater  than  1 .4  units,  so  shocking  is  not  responsible  for  the 
departure  from  Gaussian  attenuation.  A  possible  explanation  for  the  deviation 
may  be  that  the  noise  (4.4)  has  developed  strong  correlations.  This  is  a  subject 
of  future  investigations. 


Fig  4.1  Attenuation  of  a  signal  with  frequency  50  as  a  function  of  distance  due  to  the  noise  (4.4) 
with  amplitude  uo  =  0.01 .  The  points  are  the  result  of  numerical  integration  with  20  time 
steps  per  cycle  of  the  signal.  The  curve  corresponds  to  theory. 


B.  RESOLUTION  OF  DETECTION  AND  GAUSSIAN-TO-BESSEL 
TRANSITION 


The  result  obtained  above  depends  on  the  resolution  of  the  detection.  If 
the  averaging  time  is  less  than  the  period  of  the  lowest  frequencies,  poor 
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resolution  will  result  in  fundamentally  different  behavior.  Rather  than  the 
monotonic  decrease  on  the  amplitude,  the  signal  can  exhibit  restitution  effects. 
For  discrete  noise  with  equally  spaced  frequency  components,  therefore,  poor 
resolution  restores  reversibility. 

The  resolution  of  the  detection  explains  a  question  that  naturally  arises 
when  it  is  noted  that  Rudenko’s  Gaussian  prediction  does  not  depend  upon  the 
bandwidth  of  the  noise,  but  only  upon  the  total  average  intensity.  If  the 
bandwidth  is  reduced  while  the  intensity  is  held  constant,  the  amplitude  of  a 
signal  must  continue  to  attenuate  as  the  same  Gaussian.  In  the  limit  when  the 
bandwidth  is  zero,  however,  the  amplitude  of  a  signal  must  have  a  Bessel- 
function  dependence  (Ch.  II,  Sec.  E  and  Ch.  Ill,  Sec.  C-2).  How  is  the  transition 
from  Gaussian  to  Bessel  dependence  made?  In  principle,  as  long  as  the 
bandwidth  is  not  exactly  zero,  the  signal  will  attenuate  as  the  Gaussian,  although 
this  will  only  be  apparent  for  averaging  times  that  are  substantially  greater  than 
the  inverse  of  the  bandwidth.  In  practice,  the  resolution  of  the  detection  will 
eventually  be  insufficient  as  the  bandwidth  is  reduced.  The  signal  will  then 
effectively  be  in  the  presence  of  a  pure  tone,  and  the  Bessel  dependence  of  the 
signal’s  amplitude  will  be  observed.  This  is  a  striking  situation  in  which  the  time 
reversibility  of  a  system  depends  upon  the  observer.  That  is,  the  fundamental 
property  of  time  reversibility  can  only  be  established  if  the  overall  system  of  the 
waves  and  the  observer  is  considered. 

Figure  4.2  shows  the  attenuation  of  a  signal  with  frequency  50  due  to  the 

noise 
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for  two  different  resolutions.  The  band  of  frequencies  is  between  1.00  and  1.25, 
with  N  =  50  components  and  amplitude  Uo  =  0.01 .  The  results  were  obtained 
with  20  time  steps  per  cycle  of  the  signal. 
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Fig  4.2  Attenuation  of  a  signal  with  frequency  50  as  a  function  of  distance  due  to  the  noise  (4.7) 
with  amplitude  Uo  =  0.01 .  The  circles  are  the  result  of  numerical  integration  with  20  time 
steps  per  cycle  of  the  signal  and  averaging  over  a  time  of  400  (one  period  of  the  step 
frequency  0.0025).  The  squares  correspond  to  averaging  over  a  time  of  1  (one  period  of 
the  lowest  frequency  1  of  the  band).  The  curve  corresponds  to  the  Gaussian  attenuation 
theory. 

The  results  of  high  resolution  are  represented  by  circles.  In  this  case,  the 
averaging  time  is  400.  This  time  equals  a  period  of  the  step  frequency  0.0025, 
which  is  the  smallest  frequency  resulting  from  nonlinear  interactions.  The  curve 
represents  the  Gaussian  attenuation  theory.  There  is  good  agreement  because 
the  sidebands,  which  result  from  different  components  of  the  noise  interacting 
with  both  themselves  and  the  signal,  are  all  averaged  out  (refer  to  Sec.  B). 

The  squares  correspond  to  averaging  over  a  time  equal  to  one  period  of 
the  lowest  frequency  of  the  band.  This  is  the  low  resolution  limit.  The  restitution 
effects  of  the  signal  are  apparent.  The  location  of  the  zero  depends  on  the 
effective  intensity  of  the  noise  when  there  is  a  lack  of  resolution  of  the 
components  of  the  noise,  and  thus  may  depend  on  the  phases  of  the 
components.  Future  investigations  will  address  this  issue. 
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C.  DOWNSTREAM  INJECTION  OF  A  SIGNAL 


As  a  test  of  the  breakdown  of  translational  invariance  we  consider  the 
evolution  of  a  signal  that  is  injected  downstream  relative  to  the  noise  at  the 
origin.  The  signal  is  expected  to  attenuate  monotonically  as  a  Gaussian.  This  is 
suggested  by  the  independence  on  the  spectral  shape  of  the  noise  in  Rudenko’s 
solution  where  the  attenuation  coefficient  (4.3)  depends  only  on  the  total  energy 
of  the  noise.  As  the  noise  propagates  downstream,  nonlinearities  modify  the 
spectrum  but  the  total  energy  of  the  noise  remains  the  same. 

On  the  other  hand,  the  main  assumption  in  Rudenko’s  theory  is  that  the 
noise  is  normally  distributed  at  the  origin.  As  the  noise  propagates  downstream 
it  may  develop  correlations  due  to  nonlinearities.  Deviations  from  the  expected 
Gaussian  attenuation  of  a  signal  injected  downstream  may  thus  serve  as  a  probe 
of  the  degree  of  correlations  of  the  noise. 

Figure  4.3  shows  the  evolution  of  the  amplitude  of  a  signal  injected 
downstream  at  x  =  0.3  units  for  the  same  parameters  used  in  Sec.  4.1  and  Fig. 
4.1 .  We  have  chosen  the  initial  amplitude  of  the  signal  injected  at  the 
downstream  point  to  be  equal  to  the  amplitude  that  the  signal  of  Fig  4.1  injected 
at  the  origin  would  have  at  that  point.  The  solid  line  corresponds  to  the  predicted 
Gaussian  with  the  downstream  injection  point  as  the  origin.  Also  plotted  in  the 
figure  is  the  numerical  and  analytical  attenuation  of  a  signal  injected  at  the  origin. 
As  with  the  previous  results,  the  fact  that  the  numerical  noise  is  nearly  normal, 
may  be  the  origin  for  the  disagreement  between  theory  and  numerical  results  at 
distances  larger  that  1 .3  units. 
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Normalized  Amplitude 


Fig  4.3  Attenuation  of  a  signal  of  frequency  50  as  a  function  of  distance  due  to  the  noise  (4.4) 
with  amplitude  uo  =  0.01 .  The  squares  correspond  to  a  signal  injected  downstream  at  x 
0.3.  The  circles  correspond  to  the  signal  injected  at  the  origin.  The  curves  correspond 
to  theory. 
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V.  CONCLUSIONS  AND  FUTURE  WORK 


A.  CONCLUSIONS 

We  have  numerically  investigated  the  attenuation  of  a  monochromatic 
signal  in  the  presence  of  discrete  noise  in  one  dimension.  The  predicted 
Gaussian  attenuation  has  been  verified  by  the  numerical  program,  which  is 
based  on  Riemann’s  implicit  solution  of  the  exact  equation  for  the  unidirectional 
propagation  of  shockless  sound.  The  program  was  checked  for  two  cases  in 
which  analytic  solutions  exist:  propagation  of  a  pure  tone  and  the  suppression 
of  signal  by  a  pump  wave. 

We  also  presented  two  new  results.  In  the  first,  we  observed  a  transition 
from  Gaussian  to  Bessel  dependence  as  a  function  of  resolution  in  the  detection 
of  a  signal.  This  results  shows  that  the  time  reversibility  of  a  system  depends 
upon  the  observer.  That  is,  the  fundamental  property  of  time  reversibility  can 
only  be  established  if  the  overall  system  of  the  waves  and  the  observer  is 
considered. 

In  the  second  result,  we  investigated  the  evolution  of  the  amplitude  of  a 
signal  injected  downstream  from  the  noise.  We  again  observed  the  Gaussian 
attenuation.  This  result  explicitly  shows  that  the  attenuation  length  depends  on 
the  distance  the  signal  has  traveled,  thus  displaying  memory  and  breakdown  of 
translational  invariance. 


B.  FUTURE  WORK 

Throughout  this  work,  we  have  argued  that  small  but  significant 
discrepancies  between  numerical  results  and  the  analytical  predictions  may 
have  their  origin  in  the  fact  that  the  numerical  noise  is  nearly  normal.  In  fact. 
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different  sets  of  values  of  the  random  phases  yield  slightly  different  results  about 
an  average  value  that  agrees  with  theory  except  for  large  distances. 

Future  investigations  should  test  the  statistics  of  the  numerical  noise  and 
quantify  the  deviations  from  normal  distribution  at  the  origin.  Deviations  from  the 
normal  distribution  should  also  be  quantified  downstream  where  strong 
correlations  may  develop  due  to  nonlinearities.  Light  may  be  shed  onto  the 
problem  by  using  noise  with  quasiperiodic  frequencies,  where  the  spacing 
between  frequency  components  is  incommensurate.  In  principle,  nonlinearities 
cause  an  infinite  number  of  frequencies  to  be  produced  in  this  case.  Hence,  it 
may  be  that  substantially  fewer  initial  frequency  components  of  the  noise  are 
required  to  yield  an  approximate  Gaussian  attenuation  of  a  signal.  This  may  be 
important  in  applications  where  an  efficient  suppression  of  high  frequency  sound 
is  desired. 

Another  topic  of  future  investigation  should  be  stochastic 
quasimonochromatic  noise.  In  this  case,  both  the  dependence  on  the  statistics 
of  the  noise  and  the  transition  from  Gaussian  to  Bessel  dependence  can  be 
probed.  Based  on  our  observation  that  this  transition  can  be  induced  by 
insufficient  resolution  of  the  detection,  we  believe  that  the  slightest  stochasticity 
of  the  noise  will  cause  the  signal  to  attenuate  as  a  Gaussian  if  the  detection  is 
accurate. 
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APPENDIX  A.  FORWARD  TIME  NUMERICAL  PROGRAM 


/*  Numerical  integration  of  Riemann*s  equation  */ 

/*  for  propagation  of  sound  in  one  dimension  */ 

/*  */ 
/★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★■A:**********/ 

/*  This  program  simulates  the  propagation  of  a  CW  in  the  presence 
of  one  dimensional  noise  with  a  flat  distribution.  The 
signal's  amplitude  is  calculated  as  a  function  of  distance. 

The  noise  is  discrete  and  the  number  of  components  N  can  be 
made  to  vary  from  one  to  several.  Regardless  of  the  number 
of  components,  the  total  energy  of  the  noise  is  kept  constant 
by  dividing  the  peak  amplitude  of  each  noise  component  by  the 
square  root  of  N. 

LAST  UPDATE  :  21  OCTOBER  1996  */ 

/*  PROCESSOR  DIRECTIVES 

#include  "math.h" 

#include  "stdio.h" 

#include  "stdlib . h" 

/*  MACRO  DEFINITIONS  */ 

#define  pi  3.1415926535898 
ttdefine  SIZE  100002 


main ( ) 

{ 

int 

N, 

i/ j/k,l; 

/*  number  of  noise  components  */ 

double 

U[SIZE]  , 
u[SIZE] , 
T[SIZE] , 
tau[SIZE] 
r[102], 
R[203]  , 
Rl[203] , 

/*  dimensionless  particle  velocity  arrays  */ 

/*  dimensionless  time  arrays  */ 

/ 

/*  random  seed  */ 

uO, 

us, 

fo. 

/^dimensionless  peak  velocity  of  one  of 
components  */ 

/^dimensionless  peak  velocity  of  one  of 
components  */ 
/^frequency  of  signal  ^/ 

noise 

noise 

X,  dx, 
t,  dt. 

/^distance  and  distance  step*/ 

/*time  and  time  step*/ 

max, 

seed. 

/*maximum  number  generated  from  random 
/*seed  for  random  number  generation*/ 

function*/ 

dl,d2, 
al, a2, 
kl,k2,k3. 

/*variables  used  for  the  LOCK-IN*/ 
k4; 

/*  Reading  the  data  form  the  input  file  */ 
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scanf ("%lf 


%lf  %lf  %lf  %d\n", &uO, &US , &f 0, &seed, &N) ; 


/*  Variables  Initialization  */ 

max=pow (2 .0,31.0)-1.0; 
dl=d2=t=x=0 . 0 ; 

uO=uO*sqrt ( 1 . 0/N) ;  /*thus,  the  noise  energy  is  the  input 

u0^2  */ 

dt=0.001; 
dx=0 . 1 ; 


for(k=l;  k<30;  k++) {  /^Initialize  R  matrices  to  zero  */ 

R[k]=0.0; 

} 

for (1=0;  1<20;  1++) { 
srandom(seed)  ; 

for(i=0;  i<N;  i++) {  /^Random  phase  of  noise  components  */ 

r [i] = (2 . 0*pi*random( ) /max) ; 

} 

for{k=l;  k<20;  k++) { 

for(i=0;  i<SIZE;  i++) { 

U[i]=T[i]=0.0; 

} 

for(i=0;  i<N;  i++) {  /*Flat  noise  in  a  band  between 

1  and  2  */ 

for(j=0;  j<SIZE;  { 

T[j]=j*dt; 

u  [  j  ]  =uO*sin  (2 . 0*pi*  (1.0+  (double)  (i) 

/ (double) (N) ) *T[ j] +r [i] ) ; 

U[j]=U[j]+u[j] ; 

1 

) 


/*Propagate  up  to  dist  of  k{inax)*dx*/ 

/*Time  and  velocity  arrays 
initialization  */ 


for(i=0;  i<SIZE;  i++) {  /*Boundary  conditions*/ 
U[i]=U[i]+us*sin(2.0*pi*f0*T[i] )  ; 

} 


for(i=0;  i<SIZE;  i++) {  /*Rieraann's  propagation  time*/ 
tau [i] =T [i] -x; 
tau [ i ] =tau [ i ] -x/ { 1 . 0+U [ i ] ) ; 

} 


for(i=0;  i<SIZE;  i++) {  /*Extracts  the  signal  component  only*/ 
U[i]=sin(2.0*pi*f0*tau[i]  )  ; 

) 


for(i=0;  i<(SIZE-2);  i++) { /*In-phase  and  quadrature  components*/ 

al=U[i];  a2=U[i+l]; 
t=i*dt; 
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kl=cos (2 . 0*pi*f 0*t) ; 
k2=cos {2 . 0*pi*f 0* (t+dt ) ) ; 
k3=sin(2.0*pi*f0*t) ; 
k4=sin (2 . 0*pi*f 0* (t+dt) )  ; 

dl=dl+( ( (al*kl)+(a2*k2) ) * ( tau [i+1] -tau [i] ) /2.0) ; 
d2=d2+( ( (al*k3)+(a2*k4) ) * (tau [i+1] -tau [i] ) /2 . 0) ; 

} 

R1 [k] =2 . 0*f 0*sqrt (dl*dl+d2*d2 ) / (5000.0) ; 
x=k*dx;  dl=d2=t=0 . 0; 

} 


seed=seed+l . 0 ; 
for(k=l;  k<20;  k++) { 

R[k]=R[k]+Rl[k]  ; 

} 

x=0 . 0; 

} 

x=0 . 0; 

for(k=l;  k<20;  k++) {  /*  Print  the  result  */ 

printf  (**%f  %.12f\n",x,R[kl  /20.0)  ; 

x=k*dx; 

} 
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APPENDIX  B.  FORWARD  TIME  NUMERICAL  PROGRAM:  DOWNSTREAM 

INJECTION 


/★****★**★***************★*****★*★★**★★**★★★****>*•★•*•****■*■**•★■*•★•*•★*★**★•*•* 
/*  Numerical  integration  of  Riemann*s  equation  */ 

/*  for  propagation  of  sound  in  one  dimension  */ 

/*  */ 

/*  This  program  simulates  the  propagation  of  a  CW  in  the  presence 
of  one  dimensional  noise  with  a  flat  distribution.  In  particular, 
this  is  designed  to  see  how  the  absorption  of  signal  put  into  the 
developed  noise  works. 

The  signal's  amplitude  is  calculated  as  a  function  of  distance. 
The  noise  is  discrete  and  the  number  of  components  N  can  be 
made  to  vary  from  one  to  several.  Regardless  of  the  number 
of  components,  the  total  energy  of  the  noise  is  kept  constant 
by  dividing  the  peak  amplitude  of  each  noise  component  by  the 
square  root  of  N. 

LAST  UPDATE  :  21  OCTOBER  1996  */ 

/*  PROCESSOR  DIRECTIVES  */ 

#include  "math . h" 

#include  "stdio . h" 

#include  "stdlib.h" 

/*  MACRO  DEFINITIONS  */ 

#define  pi  3.14159265 
ttdefine  SIZE  100002 


main ( ) 

{ 

int 

N, 

i, j,k,l; 

/*  number  of  noise  components  */ 

double 

U[SIZE] , 
u[SIZE]  , 
T[SIZE]  , 
tau[SIZE] , 
r[102], 
R[203]  , 
Rl[203] , 

/*  dimensionless  particle  velocity  arrays  */ 

/*  dimensionless  time  arrays  */ 

/*  random  seed  */ 

uO, 

us, 

fo. 

/^dimensionless  peak  velocity  of  one  of  noise 
components  */ 

/^dimensionless  peak  velocity  of  one  of  noise 
components  */ 

/^frequency  of  signal  */ 

X,  dx, 
t,  dt. 

/*distance  and  distance  step*/ 

/*time  and  time  step*/ 

max, 

seed. 

/*maximum  number  generated  from  random  function*/ 
/*seed  for  random  number  generation*/ 

dl,d2, 
al, a2. 

/*variables  used  for  the  LOCK-IN*/ 
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kl,k2,k3,k4; 

/*  Reading  the  data  form  the  input  file  */ 

scanf("%lf  %lf  %lf  %lf  %d\n”,  &u0,  &US;. &f  0,  &seed,  &N)  ; 

/*  Variables  Initialization  */ 

max=pow (2 .0,31.0)-1.0; 
dl=d2=t=x=0, 0; 

u0=u0*sqrt (1 . 0/N) ;  /*thus,  the  noise  energy  is  the  input 

u0^2  */ 

dt=0.001; 
dx=0 . 1; 

us=us*exp (” (2 . 0*pi*50 . 0*u0*0 , 3) * (2 . 0*pi*50 . 0*u0*0 . 3) *N/4 . 0 ) ; 

for(k=l;  k<30;  k++) {  /^Initialize  R  matrices  to  zero*/ 

R[k]=0.0; 

} 

for (1=0;  1<20;  1++) { 
s  random  (seed)  ; 

for(i=0;  i<N;  i++) {  /*Random  phase  of  noise  components*/ 

r [i] = (2 . 0*pi*random( ) /max) ; 

} 

for(k=3;  k<20;  k++) {  /*Propagate  up  to  dist  of  k(max)*dx*/ 


/*  Propagate  the  noise  componets  up  to  0.3m  first  */ 


x=0 . 3; 

for(i=0;  i<SIZE;  i++) {  /*Time  and  velocity  arrays 

initialization  */ 

U[i]=T[i]=0. 0; 

) 

for{i=0;  i<N;  i++) {  /*Flat  noise  in  a  band  between 

1  and  2*/ 

for(j=0;  j<SIZE;  j++)  { 

T[j]=j*dt; 

u  [  j  ]  =u0*cos  (2 . 0*pi*  (1.0+  (double)  (i)  / 

(double) (N) ) *T [ j ] +r [i] ) ; 

U[j]=U[j]+u[j]  ; 

} 

} 


for(i=0;  i<SIZE;  i++) {  /*Riemann's  propagation  time  for 

noise*/ 

T[i]=T[i]-x/ (1.0+U[i] )  ; 

} 


/*  Add  signal  to  developed  noise  and  propagate  both  of  them  */ 

/*  up  to  k(max)*dx  */ 
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x=k*dx; 

for(i=0;  i<SIZE;  i++) {  /*Inject  the  signal  */ 
U[i]=U[i]+us*cos (2.0*pi*fO*T[i] ) ; 

) 


for(i=0;  i<SIZE;  i++) {  /*Riemann's  propagation  time  for  signal 

and  noise*/ 

tau [i]=T[i]; 

tau[i]=tau [i] - (x-0. 3) / { 1 . 0+U [i] ) ; 

) 


for{i=0;  i<SIZE;  i++) {  /*Extracts  the  signal  component  only*/ 
U[i]=us*cos (2 . 0*pi*f0*T [i] ) ; 

) 

for(i=0;  i<SIZE-2;  i++) {  /*In-phase  and  quadrature  components*/ 

al=U[i];  a2=U[i+l]; 
t=i*dt; 

kl=cos (2 . 0*pi*f0*tau[i] ) ; 
k2=cos (2 . 0*pi*f0* (t+dt) ) ; 
k3=sin(2.0*pi*f0*tau[i] ) ; 
k4=sin (2 . 0*pi*f0* (t+dt) ) ; 

dl=dl+(al*kl)* (tau[i+l]-tau[i] ) ; 
d2=d2+ (al*k3) * (tau [i+1] -tau [i] ) ; 


R1 [k]=2.0*f0*sqrt (dl*dl+d2*d2) / (5000.0*0.0001) ; 
x=k*dx;  dl=d2=t=0.0; 

} 

seed=seed+l . 0 ; 
for(k=3;  k<20;  k++) { 

RCk]=R[k]+Rl[k]  ; 

) 


for(k=3;  k<20;  k++) {  /*Print  the  result*/ 

x=k*dx; 

print f ("%f  % . 12f\n",x,  R[k] /20. 0)  ; 


} 


) 
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